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1. Introduction 

Disconnected quark line diagrams have non-trivial contribution to the lattice QCD study of 
nucleon structure. However, direct calculation of the disconnected diagrams is not practical as they 
require all-to-all propagators. Development of stochastic estimation for the disconnected diagrams 
made it feasible with current computing power, but it is still one of the most expensive computations 
in lattice QCD [1]. 

The simplest observable are the local bilinear operators that give the nucleon chai'ges gA,sj^ 
nucleon a term and strangeness. Those play an important role in dark matter search and study of 
neutron electric dipole moment (nEDM) [2]. In this paper we report a lattice QCD calculation of 
the disconnected quark loop contribution to isoscalar nucleon charges of light and strange quarks. 
The calculation is done with clover valence quarks on the Nf = 2 + 1 + 1 highly improved staggered 
quarks (HISQ) ensembles generated by MILC collaboration [3]. The three ensembles that we use 
are listed in Table 1 . 


ID 

a (fm) 

(MeV) 

+3 X r 

M^L 

al2m310 

0.1207(11) 

305.3(4) 

243 X 64 

4.54 

al2m220 

0.1184(10) 

216.9(2) 

32^ X 64 

4.29 

a09m310 

0.0888(08) 

312.7(6) 

323 X 96 

4.50 


Table 1: MILC HISQ lattices analyzed in this study. Lattice parameters are obtained from Ref. [3]. 


Nucleon charge of quark flavor q and gamma structure T can be extracted from the ratio of 
the three-point function and the two-point function of nucleons as follows: 




{Tv[^rCp%^ /»T»0 q 

(Tr['^2ptC2pH0]) 


( 1 . 1 ) 


Here T is the gamma structure of the local bilinear operator in the three point function, t is the 
source-sink separation, T is the distance in Euclidean time between the source and inserted operator, 
and and are the projection operators. There are two possible classes of contractions when 
we insert a bilinear operator between two nucleon states. One is to contract the bilinear operator 
with one of the three quarks of nucleons, as shown in the left hand side (l.h.s) diagram of Eig. 1, and 
the other is to contract the bilinear operator by themselves making a quark loop, as shown in the 
right hand side (r.h.s) diagram of Eig. 1. Here we report our recent calculation of the disconnected 
quark loop contributions to the nucleon charges with various improvement techniques. 


2. Disconnected Quark Loop Contribution 


The disconnected quark loop part of Eq. (1.1) is written as 


4^E^) = (LTr[M->(T,x;T,x)r]) 

X 


(Tr[,^rC^Pt(0] i:,Tr[M-‘(T,x;T,x)r]) 
(Tr[,^2ptC2pt(t)]) 


( 2 . 1 ) 


with the Dirac operator M. The calculation needs the evaluation of (1) nucleon two-point functions 
C2p‘ and (2) quark loops Tr[M^’r]. 
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Figure 1: Quark line connected (left) and disconnected (right) diagrams of the nucleon three-point function. 


2.1 Two-point Function 

Exploiting translation symmetry of lattice, we average the nucleon two-point function over 
multiple source positions. In order to reduce the computational cost for the multiple source posi¬ 
tions, we apply the all-mode-averaging (AMA) technique [4]. Using AMA, the two-point functions 
for most of the source positions are calculated by a low-precision (LP) Dirac inverter with cheap 
computational cost. The systematic error of the LP estimate is corrected by a few high-precision 
(HP) calculations of those as follows: 


1 


Nlp 


c2p‘^-p(t,to) = ^£cfp‘ 

Alp 


Lp(^;fo,x, 


1 Nup 


(t; to,xf- Clp (t; to, x} 


-,2pt / 


HP^ 


( 2 . 2 ) 


where Cpp and C^p are the two-point function calculated in LP and HP, respectively, and 
X™ are the different source positions, and App (Ahp) is the number of source positions where the 
two-point function is evaluated in LP (HP). 

In this study, we spread 60 LP and 4 HP source positions in four timeslices, and calculate LP 
estimate of the two-point function by truncating the Dirac inversion with low-accuracy stopping 
criterion, rpp ~ 10^^. Since we need to perform multiple Dirac inversions on a lattice, deflating 
the low-eigenmodes is an efficient way of reducing the computational cost. Here we use multigrid 
solver [5] to deflate the low-eigenmodes. 


2.2 Quark Loop 

The quark loops ^xTr[(H^^(T,x; T,x)r] are estimated by a stochastic method. Let us consider 
random noise vectors I?],) in color/spin/spacetime space, satisfying 


N 


Vn 


N 




A: 


Vn 


(2.3) 


Having the corresponding solutions \si) of the Dirac equation M\si) = |t];), one can stochastically 
estimate the inverse Dirac matrix as well as the quark loop by 


1 


N 




r=l 


(v^)’ 


N 


£Tr[M i(T,x;T,x)r] - £(T],-|^r|5,-)^, 


(2.4) 


r=l 
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where \x)^ is a vector whose t = x components are filled by |v), and other components are zero. 

For the stochastic estimation, we adopt the truncated solver method (TSM) [1] to reduce com¬ 
putation time. Similar to the AMA, which we use in the two-point function calculation, most of 
the calculation is done in LP, and the result is corrected by small number of HP calculation: 


Vlp 2 ^lp+Mip 

Xl (I‘^')hp- kf)Lp) {rii\ ■ 

Alp Ahp 


(2.5) 


Here |5 ,)lp and |5 ,)hp are the LP and HP solution of the Dirac equation for a given source |t],). 

For the random noise source, we tested various type of random numbers such as Z 2 , Z 2 0 /Z 2 , 
Z 4 and Gaussian. It turns out that all of them are similar in terms of the total statistical error, but 
the Gaussian random number is marginally better in our case. However, if one does not use the 
hopping parameter expansion, which will be discussed later in this section, type random noises 
will be a better choice as they cancel part of leading noises. 

The scaling form of the total statistical error of disconnected quark loop is [4] 


where /?lp/hp = Alp/Ahp, and C is determined by the difference (correlation) between |5,)hp and 
ki)LP- In other words, the error scales as a/1/Alp, and the contribution from the correction term 
is determined by /?lp/hp and accuracy of the LP estimate. In order to maximize computation 
efficiency, we use /?lp/hp = 30, and calculate |5,)lp by the multigrid solver with plp ~ 5 x 10^^. 

Required App is determined by the scaling of total eiTor in disconnected contribution to the 
nucleon charges, which is a combination of the error from the nucleon two-point function and 
that of the disconnected quark loop. App dependence of the total error can be written as a = 
CJooY^l +Atsm/Alp. Here Oco is the error when App —)• so it is determined by the precision 
of the two-point function. Xtsm is determined by the gamma structure of the operator and TSM 
parameters. On al2m310 lattice, with all the error reduction techniques given in this paper, we find 
that Xtsm ~ 40,650 and 15000 for gs, gA and gr, respectively, for light quarks. In other words, 
1 ^( 100 ) of random noise sources is enough for the scalar, but ^( 1000 ) ~ ^( 10000 ) is needed 
for the tensor. In this study, we use App = 5000(1500], 11000(4000] and 4000(1200] on al2m310, 
al2m220 and a09m310 lattices for light [strange] quarks. 

For error reduction, we use the hopping parameter expansion (HPE) [ 6 , 7]. The Dirac matrix, 
written in the form of M = ^ (1 — kD) with a hopping parameter K can be expanded as 

1 1 

—=l + Y,{KDy + {KDf—M-^, (2.7) 

2k 2k 


up to any order of k. Taking k = 2, the disconnected quark loops are given as 


Tr[M^ir] =Tr[{2Kl + 2K'^D + K^D^M-^)r] . 


( 2 . 8 ) 


Here the first term inside the trace of r.h.s can be easily calculated by hand, and the second term 
cannot contribute to the nucleon charges. Hence what we need to evaluate is only the third term, 
Tr . Since noise terms are removed, HPE allows us to reduce the statistical error. 
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One of the major sources of the statistical noise in the stochastic estimation of the discon¬ 
nected quark loops is the correlation between the spacetime/spin/color components. Hence one 
may reduce the noise by ( 1 ) dividing those components into m subspaces with proper separations, 
(2) evaluating the quark loops in those subspaces and (3) obtaining the full results by summing the 
results of those subspaces. This is the dilution technique [ 8 , 9, 10]. In dilution, the computation 
cost increases as it needs to be evaluated on all m subspaces. As a results, it is useful only if the 
gain from the noise reduction is bigger than the increase of computation cost. We do not find any 
error reduction from the dilution in spin or color components on nucleon charges. Dilution in time 
direction reduces error, but the gain is similar or smaller than the increase of computation cost. 

We average all the possible combinations of the gamma structures, and forward/backward 
propagations. We also set zero for the real or imaginary components that should be zero by the 
75 -hermiticity of clover Dirac operator, 

3. Results and Conclusion 

The nucleon charges should be calculated with the ground state nucleons. On the lattice, 
however, nucleon state is generated by an interpolating operator, and it also generates the excited 
states of the nucleon. We remove the contamination from the excited states by fitting the results to 
a fit function including terms up to one excited state: 

C^^\tf,ti) =|M)|^“0(9-0 + |j^i|, (3.1) 

where ti = 0, tf = t, and | 0 ) and | 1 ) are the ground and first excited nucleon states, respectively. 
Details of the fitting procedure is described in Ref. [2] as the two-simRR method. Extrapolation 
plots on al2m310 ensemble are given in Fig. 2. 

The extrapolated results are renormalized at 2 GeV in MS scheme by using the non-perturbative 
renormalization (NPR) only with the connected diagrams. As an approximation, disconnected dia¬ 
grams and possible mixing with other quark flavors are neglected [11]. The quark mass difference 
befween lighf and sfrange quarks is also neglected; fhe charges for fhe sfrange quarks are renormal¬ 
ized by fhose of lighf quarks. The NPR procedure is described in Ref. [2]. 

In Table 2 we presenf our preliminary resulfs of fhe connected and disconnecfed quark loop 
confribufions fo isoscalar nucleon charges of lighf quarks and nucleon charges of sfrange quarks 
on fhree differenf ensembles. The calculafion of fhe connecfed parfs is described in Ref. [2]. For 
disconnecfed we obfained slighfly smaller values (abouf 2 a) fhan fhose from aufhors of 

Ref. [12]. Disconnecfed g’^‘^ is yef fo be compared as fhey have large sfafisfical error. For g^, we 
obfained consisfenf value wifh Refs. [11, 12, 13] on a09m310 lattice, and slighfly smaller value on 
al2 lattices. However, nofe fhaf confinuum and chiral exfrapolafions are needed fo give a meaning¬ 
ful comparison. 

In fhis sfudy, we calculate fhe disconnecfed quark loop confribufion fo fhe nucleon charges. 
Those are imporfanf observables for fhe neufron elecfric dipole momenf (nEDM) sfudy fhrough fhe 
quark EDM, and dark maffer search. The fechniques of disconnecfed quark loops presenfed in fhis 
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x/a 



3 4 5 6 7 8 9 

x/a 

Figure 2: Disconnected contribution to the nucleon charges for different source-sink separation fsep and 
operator insertion x in lattice unit on al2m310 ensemble. Those are extrapolated to the limit of fsep ^ T ^ 0, 
plotted as a black solid line with gray error band. 


paper can be applied to the study of hadronic form factors, chromo EDM operators and transverse 
momentum dependent distribution function (TMD) observables [14]. 
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Charge 

al2m310 

al2m220 

a09m310 

Conn. 

Diseonn. 

Conn. 

Diseonn. 

Conn. 

Diseonn. 

sf" 

5.22(23) 

1.65(14) 

6.16(33) 

1.49(36) 

6.01(25) 

1.40(14) 

sU 

0.581(25) 

-0.116(12) 

0.587(31) 

-0.170(41) 

0.628(22) 

-0.103(16) 


0.619(26) 

-0.024(05) 

0.625(31) 

-0.005(11) 

0.622(25) 

-0.011(05) 


0 

0.584(48) 

0 

0.639(81) 

0 

0.440(44) 


0 

-0.036(05) 

0 

-0.038(10) 

0 

-0.022(07) 

St 

0 

-0.0027(24) 

0 

-0.0009(32) 

0 

0.0038(34) 


Table 2: Preliminary results of the connected and disconnected quark loop contributions to isoscalar nucleon 
charges of light quarks and nucleon charges of strange quarks on three different ensembles. Those are 
renormalized at 2GeV in MS scheme, using NPR only with the connected part. 
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